Before we jump into any statistical analyses, just want to get an idea of how things play out by treatment and temperature .
| tempxtrt | treatment | temperature | meanchangeTDN_uM | sdchangeTDN_uM | n |
|---|---|---|---|---|---|
| 8Control | C | 8 | 7.204 | 2.053 | 6 |
| 8N+P | NP | 8 | 48.590 | 3.815 | 6 |
| 8Nitrogen | N | 8 | 14.838 | 3.488 | 6 |
| 8Phosphorus | P | 8 | 9.863 | 0.542 | 5 |
| 12Control | C | 12 | 10.255 | 0.826 | 6 |
| 12N+P | NP | 12 | 58.205 | 1.072 | 6 |
| 12Nitrogen | N | 12 | 23.563 | 5.616 | 6 |
| 12Phosphorus | P | 12 | 10.051 | 0.805 | 6 |
| 16Control | C | 16 | 8.045 | 1.400 | 6 |
| 16N+P | NP | 16 | 54.058 | 1.616 | 6 |
| 16Nitrogen | N | 16 | 43.418 | 3.416 | 6 |
| 16Phosphorus | P | 16 | 8.834 | 1.614 | 6 |
Here we have some summary statistics for TDN uptake. First we are looking at change in TDN from start to end of experiment by TEMPERTURE and by TREATMENT separately. I’ve also plotted boxplots with individual points.
## -------------------------------------------------------------------------
## changeTDN_uM ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean 2.1e+01 2.6e+01 2.9e+01
## median 1.1e+01 1.4e+01 2.5e+01
## sd 1.7e+01 2.0e+01 2.1e+01
## IQR 2.2e+01 2.7e+01 4.0e+01
## n 2e+01 2e+01 2e+01
## np 32.4% 33.8% 33.8%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 2.5284, df = 2, p-value = 0.2825
## -------------------------------------------------------------------------
## changeTDN_uM ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N P NP
## mean 8.5e+00 2.7e+01 9.6e+00 5.4e+01
## median 8.8e+00 2.1e+01 9.9e+00 5.4e+01
## sd 1.9e+00 1.3e+01 1.2e+00 4.7e+00
## IQR 3.6e+00 2.3e+01 1.5e+00 5.7e+00
## n 2e+01 2e+01 2e+01 2e+01
## np 25.4% 25.4% 23.9% 25.4%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 58.424, df = 3, p-value = 1.276e-12
# Change stripchart colors by groups
ggplot(expchem_trim, aes(x=temperature, y=changeTDN_uM, color=treatment)) +
geom_boxplot(position=position_dodge(0.8))+
geom_jitter(position=position_jitter(0.2))+
scale_color_manual(values=c("#40A339", "#3166A9", "#DE0019", "#82388E"))+
theme_classic2(base_size=16)
ggplot(expchem_trim, aes(x=temperature, y=changeTDN_uM, color=treatment)) +
geom_boxplot(position=position_dodge(0.8))+
labs(x="Temperature Treatment",y="Nitrogen Uptake (uM)")+
scale_color_manual(values=c("#40A339", "#3166A9", "#DE0019", "#82388E"))+
theme_bw(base_size=16)
trtcolors <- c("#40A339", "#3166A9", "#DE0019", "#82388E") #treatment colors
TDN_uptakerate <- ggplot(expchem_changesummary, aes(x=treatment, y=meanTDNuptake_dailyrate, fill=treatment)) + theme_bw()+
geom_bar(color="black",stat="identity", position="dodge") +
scale_fill_manual(values=trtcolors,
name="Nutrient\namendment",
labels=c("Control", "Nitrogen", "Phosphorus", "Both (N & P)"))+
geom_errorbar(aes(ymin=meanTDNuptake_dailyrate-sdTDNuptake_dailyrate, ymax=meanTDNuptake_dailyrate+sdTDNuptake_dailyrate),
width=.2,
position=position_dodge(0.9))+
coord_cartesian(ylim=c(0,8))+
scale_y_continuous(breaks=seq(0,8,1))+
theme(axis.text.x = element_text(colour = "black", size = 40), axis.title.y=element_text(size=36, vjust=1.2, face="bold"),
axis.text.y = element_text(colour = "black", size = 30), axis.title.x=element_text(size=36, vjust=-0.1, face="bold"),
legend.text = element_text(colour="black", size = 30), legend.title=element_text(size=36, face="bold"),
plot.title=element_text(size=40, face="bold"),
panel.border = element_rect(size=1.0, fill=NA),
legend.key.size = unit(3, 'lines'),
panel.grid.minor = element_blank(),
panel.grid.major = element_blank(),
strip.background =element_rect(fill="white",linetype = 0),
strip.text = element_text(colour = 'black', size=28))+
facet_wrap( ~ temperature , ncol=3, scales="free")+
# geom_hline(yintercept=66) +
# annotate("text",1,66,vjust=-1,label=" Starting [TDN] in N & N+P treatments")+
# geom_hline(yintercept=18, color="black", linetype="dashed")+
# annotate("text",0.8,18,vjust=-1,label="Background [TDN]")+
labs(x="",y=" [uM N/day]", title="Dissolved inorganic N uptake rate\n")
TDN_uptakerate
# ggsave("TDN_uptakerate.png", width=16, height=9.5,units="in")
trtcolors <- c("#40A339", "#3166A9", "#DE0019", "#82388E") #treatment colors
TDN_uptakerate_boxplot <- ggplot(expchem_trim, aes(x=temperature, y=TDNuptake_dailyrate, fill=treatment)) +
geom_boxplot()+
stat_boxplot(geom ='errorbar') +
scale_fill_manual(values=trtcolors,
name="Nutrient\nTreatment",
labels=c("Control", "Nitrogen", "Phosphorus", "N&P"))+
coord_cartesian(ylim=c(0,8))+
scale_y_continuous(breaks=seq(0, 8, 1))+
theme(axis.text.x = element_text(colour = "black", size = 40), axis.title.y=element_text(size=36, vjust=1.2, face="bold"),
axis.text.y = element_text(colour = "black", size = 30), axis.title.x=element_text(size=36, vjust=-0.1, face="bold"),
legend.text = element_text(colour="black", size = 30), legend.title=element_text(size=36, face="bold"),
legend.key.size = unit(3, 'lines'),
plot.title=element_text(size=40,face="bold"),
panel.border = element_rect(size=1.0, fill=NA),
panel.grid.minor = element_blank(),
panel.grid.major = element_blank())+
labs(x="Incubation temperature",y="[uM/day]", title="Dissolved inorganic N uptake rate")
TDN_uptakerate_boxplot
# ggsave("TDN_uptakerate_boxplot.png", width=16, height=9.5,units="in")
#What percentage of N was taken up
ggplot(expchem_trim, aes(x=treatment, y=percTDNuptake, color=treatment)) +
geom_boxplot(position=position_dodge(0.8))+
geom_jitter(position=position_jitter(0.2))+
scale_color_manual(values=c("#40A339", "#3166A9", "#DE0019", "#82388E"))+
theme_classic2(base_size=16)+
facet_wrap( ~ temperature, ncol=3, scales="free")+
labs(x="Nutrient Treatment",y="Percent of Initial TDN assimilated")
#Bar plot of % N taken up
expcolors <- c("#40A339", "#3166A9", "#DE0019", "#82388E")
percNuptake <- ggplot(expchem_changesummary, aes(x=treatment, y=meanpercTDNuptake, fill=treatment)) + theme_bw()+
geom_bar(color="black",stat="identity", position="dodge") +
scale_fill_manual(values=expcolors,
name="Nutrient\nTreatment\n",
labels=c("Control\n", "Nitrogen\n", "Phosphorus\n", "Both (N & P)\n"))+
geom_errorbar(aes(ymin=meanpercTDNuptake-sdpercTDNuptake, ymax=meanpercTDNuptake+sdpercTDNuptake),
width=.2,
position=position_dodge(0.9))+
coord_cartesian(ylim=c(0,100))+
scale_y_continuous(breaks=seq(0,100,10))+
theme(axis.text.x = element_text(colour = "black", size = 28), axis.title.y=element_text(size=30, vjust=1.2, face="bold"),
axis.text.y = element_text(colour = "black", size = 28), axis.title.x=element_text(size=30, vjust=-0.1, face="bold"),
legend.text = element_text(colour="black", size = 28), legend.title=element_text(size=30, face="bold"),
plot.title=element_text(size=30, face="bold"),
panel.border = element_rect(size=1.0, fill=NA),
panel.grid.minor = element_blank(),
panel.grid.major = element_blank(),
strip.background =element_rect(fill="white"),
strip.text = element_text(colour = 'black', size=28))+
facet_wrap( ~ temperature, ncol=3, scales="free")+
# geom_hline(yintercept=66) +
# annotate("text",1,66,vjust=-1,label=" Starting [TDN] in N & N+P treatments")+
# geom_hline(yintercept=18, color="black", linetype="dashed")+
# annotate("text",0.8,18,vjust=-1,label="Background [TDN]")+
labs(x="Temperature Treatment",y="% Assimilated", title="Percentage of nitrogen assimilated from experimental starting conditions")
percNuptake
ggsave("percNuptake.png", width=20, height=12,units="in")
From our preliminary glance at the data, it appears there is an interaction between the temperature of the experiment and the nutrient treatment. In the N-only treatment, TDN uptake increased with temperature. Interesting, in the N+P nutrient treatment, almost all of the TDN was consumed in the experiment. The following two plots are interaction plots of by nutrient treatment and temperature treatment, respecitvely.
## Anova Table (Type III tests)
##
## Response: changeTDN_uM
## Sum Sq Df F value Pr(>F)
## (Intercept) 311.4 1 43.721 1.222e-08 ***
## temperature 29.8 2 2.093 0.1324
## treatment 6553.0 3 306.723 < 2.2e-16 ***
## temperature:treatment 1980.6 6 46.354 < 2.2e-16 ***
## Residuals 420.2 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Anova Table (Type III tests)
##
## Response: percTDNuptake
## Sum Sq Df F value Pr(>F)
## (Intercept) 9517.6 1 231.524 < 2.2e-16 ***
## temperature 911.2 2 11.083 8.193e-05 ***
## treatment 8349.3 3 67.702 < 2.2e-16 ***
## temperature:treatment 5107.9 6 20.709 6.879e-13 ***
## Residuals 2425.4 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), we see there is a evidence of a significant interaction effect (F=46.453, p<0.0001). The test for the main effect of nutrient treatment is significant (p<0.0001) but the test for the main effect of temperature alone is not significant (p=0.13). After fitting an ANOVA model it is important to always check the relevant model assumptions. This includes making QQ-plots and residual plots. Will use Levene’s test of equality of variances too.
##
## Call:
## lm(formula = changeTDN_uM ~ treatment * temperature, data = expchem_trim)
##
## Residuals:
## Min 1Q Median 3Q Max
## -6.898 -1.163 -0.040 1.105 7.392
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.2037 1.0895 6.612 1.22e-08 ***
## treatmentN 7.6347 1.5407 4.955 6.41e-06 ***
## treatmentP 2.6593 1.6159 1.646 0.10514
## treatmentNP 41.3863 1.5407 26.862 < 2e-16 ***
## temperature12 3.0517 1.5407 1.981 0.05230 .
## temperature16 0.8417 1.5407 0.546 0.58693
## treatmentN:temperature12 5.6733 2.1789 2.604 0.01165 *
## treatmentP:temperature12 -2.8640 2.2327 -1.283 0.20460
## treatmentNP:temperature12 6.5633 2.1789 3.012 0.00382 **
## treatmentN:temperature16 27.7383 2.1789 12.730 < 2e-16 ***
## treatmentP:temperature16 -1.8707 2.2327 -0.838 0.40550
## treatmentNP:temperature16 4.6267 2.1789 2.123 0.03793 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.669 on 59 degrees of freedom
## Multiple R-squared: 0.9845, Adjusted R-squared: 0.9816
## F-statistic: 340.5 on 11 and 59 DF, p-value: < 2.2e-16
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 1.8448 0.06639 .
## 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Error in library(interplot): there is no package called 'interplot'
Neither plot indicates a significant violation of normality assumptions. We’re good here.
Pairwise comparsions for each level of temperature
## $emmeans
## temperature = 8:
## treatment emmean SE df lower.CL upper.CL
## C 7.203667 1.089457 59 5.023668 9.383665
## N 14.838333 1.089457 59 12.658335 17.018332
## P 9.863000 1.193440 59 7.474931 12.251069
## NP 48.590000 1.089457 59 46.410002 50.769998
##
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C 10.255333 1.089457 59 8.075335 12.435332
## N 23.563333 1.089457 59 21.383335 25.743332
## P 10.050667 1.089457 59 7.870668 12.230665
## NP 58.205000 1.089457 59 56.025002 60.384998
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C 8.045333 1.089457 59 5.865335 10.225332
## N 43.418333 1.089457 59 41.238335 45.598332
## P 8.834000 1.089457 59 6.654002 11.013998
## NP 54.058333 1.089457 59 51.878335 56.238332
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 8:
## contrast estimate SE df t.ratio p.value
## C - N -7.6346667 1.540725 59 -4.955 <.0001
## C - P -2.6593333 1.615926 59 -1.646 0.3615
## C - NP -41.3863333 1.540725 59 -26.862 <.0001
## N - P 4.9753333 1.615926 59 3.079 0.0162
## N - NP -33.7516667 1.540725 59 -21.906 <.0001
## P - NP -38.7270000 1.615926 59 -23.966 <.0001
##
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N -13.3080000 1.540725 59 -8.637 <.0001
## C - P 0.2046667 1.540725 59 0.133 0.9992
## C - NP -47.9496667 1.540725 59 -31.121 <.0001
## N - P 13.5126667 1.540725 59 8.770 <.0001
## N - NP -34.6416667 1.540725 59 -22.484 <.0001
## P - NP -48.1543333 1.540725 59 -31.254 <.0001
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N -35.3730000 1.540725 59 -22.959 <.0001
## C - P -0.7886667 1.540725 59 -0.512 0.9559
## C - NP -46.0130000 1.540725 59 -29.865 <.0001
## N - P 34.5843333 1.540725 59 22.447 <.0001
## N - NP -10.6400000 1.540725 59 -6.906 <.0001
## P - NP -45.2243333 1.540725 59 -29.353 <.0001
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Pairwise comparsions for each level of nutrient treatment
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 7.203667 1.089457 59 5.023668 9.383665
## 12 10.255333 1.089457 59 8.075335 12.435332
## 16 8.045333 1.089457 59 5.865335 10.225332
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 14.838333 1.089457 59 12.658335 17.018332
## 12 23.563333 1.089457 59 21.383335 25.743332
## 16 43.418333 1.089457 59 41.238335 45.598332
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 9.863000 1.193440 59 7.474931 12.251069
## 12 10.050667 1.089457 59 7.870668 12.230665
## 16 8.834000 1.089457 59 6.654002 11.013998
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 48.590000 1.089457 59 46.410002 50.769998
## 12 58.205000 1.089457 59 56.025002 60.384998
## 16 54.058333 1.089457 59 51.878335 56.238332
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -3.0516667 1.540725 59 -1.981 0.1259
## 8 - 16 -0.8416667 1.540725 59 -0.546 0.8488
## 12 - 16 2.2100000 1.540725 59 1.434 0.3301
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -8.7250000 1.540725 59 -5.663 <.0001
## 8 - 16 -28.5800000 1.540725 59 -18.550 <.0001
## 12 - 16 -19.8550000 1.540725 59 -12.887 <.0001
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -0.1876667 1.615926 59 -0.116 0.9926
## 8 - 16 1.0290000 1.615926 59 0.637 0.8005
## 12 - 16 1.2166667 1.540725 59 0.790 0.7108
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -9.6150000 1.540725 59 -6.241 <.0001
## 8 - 16 -5.4683333 1.540725 59 -3.549 0.0022
## 12 - 16 4.1466667 1.540725 59 2.691 0.0247
##
## P value adjustment: tukey method for comparing a family of 3 estimates
## Error in UseMethod("cld"): no applicable method for 'cld' applied to an object of class "c('emm_list', 'list')"
## Error in UseMethod("cld"): no applicable method for 'cld' applied to an object of class "c('emm_list', 'list')"
Here we have some summary statistics for NH4 at end of the experiment.
## -------------------------------------------------------------------------
## NH4_mgL ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean 1.0e-01 1.1e-01 4.5e-02
## median 5.0e-02 6.0e-02 4.0e-02
## sd 1.1e-01 1.2e-01 2.0e-02
## IQR 1.3e-01 1.4e-01 2.3e-02
## n 2e+01 2e+01 2e+01
## np 32.4% 33.8% 33.8%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 2.8018, df = 2, p-value = 0.2464
## -------------------------------------------------------------------------
## NH4_mgL ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N NP P
## mean 5.0e-02 1.9e-01 5.2e-02 5.5e-02
## median 4.0e-02 2.3e-01 4.5e-02 3.0e-02
## sd 3.4e-02 1.0e-01 2.0e-02 1.1e-01
## IQR 2.0e-02 1.9e-01 3.0e-02 2.0e-02
## n 2e+01 2e+01 2e+01 2e+01
## np 25.4% 25.4% 25.4% 23.9%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 31.577, df = 3, p-value = 6.425e-07
# Change stripchart colors by groups
ggplot(expchem, aes(x=temperature, y=NH4_mgL, color=treatment)) +
geom_boxplot(position=position_dodge(0.8))+
labs(x="Temperature Treatment",y="Ammonia conc. at end (mg/L)")+
scale_color_manual(values=c("#40A339", "#3166A9", "#DE0019", "#82388E"))+
theme_bw(base_size=16)
ggsave("NH4_end_bytemp.png", width=10, height=6,units="in")
ggplot(expchem, aes(x=treatment, y=NH4_mgL, color=temperature)) +
geom_boxplot(position=position_dodge(0.8))+
labs(x="Temperature Treatment",y="Ammonia conc. at end (mg/L)")+
theme_bw(base_size=16)
ggsave("NH4_end_bytrt.png", width=10, height=6,units="in")
There is an effect of temperature alone and there is an effect of nutrient treatment alone and an interaction between temp*nutrients.
## Anova Table (Type III tests)
##
## Response: NH4_mgL
## Sum Sq Df F value Pr(>F)
## temperature 0.04890 3 4.4939 0.006598 **
## treatment 0.27515 3 25.2860 1.208e-10 ***
## temperature:treatment 0.13111 6 6.0244 5.818e-05 ***
## Residuals 0.21400 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
After fitting an ANOVA model it is important to always check the relevant model assumptions. This includes making QQ-plots and residual plots. Will use Levene’s test of equality of variances too.
##
## Call:
## lm(formula = NH4_mgL ~ temperature * treatment - 1, data = expchem)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.08500 -0.01833 -0.00333 0.01333 0.38500
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## temperature8 0.045000 0.024587 1.830 0.07227 .
## temperature12 0.070000 0.024587 2.847 0.00606 **
## temperature16 0.035000 0.024587 1.424 0.15986
## treatmentN 0.241667 0.034772 6.950 3.29e-09 ***
## treatmentNP -0.001667 0.034772 -0.048 0.96193
## treatmentP -0.021000 0.036469 -0.576 0.56692
## temperature12:treatmentN -0.085000 0.049174 -1.729 0.08912 .
## temperature16:treatmentN -0.218333 0.049174 -4.440 4.01e-05 ***
## temperature12:treatmentNP -0.011667 0.049174 -0.237 0.81328
## temperature16:treatmentNP 0.023333 0.049174 0.475 0.63689
## temperature12:treatmentP 0.056000 0.050389 1.111 0.27092
## temperature16:treatmentP 0.017667 0.050389 0.351 0.72713
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.06023 on 59 degrees of freedom
## Multiple R-squared: 0.8231, Adjusted R-squared: 0.7871
## F-statistic: 22.87 on 12 and 59 DF, p-value: < 2.2e-16
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 0.7057 0.7281
## 59
We don’t exactly meet the assumptions of normality.
##
## Call:
## lm(formula = logNH4 ~ temperature * treatment - 1, data = expchem)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.46007 -0.13249 -0.00846 0.12301 1.04024
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## temperature8 -1.39613 0.09382 -14.880 < 2e-16 ***
## temperature12 -1.23594 0.09382 -13.173 < 2e-16 ***
## temperature16 -1.46508 0.09382 -15.615 < 2e-16 ***
## treatmentN 0.85174 0.13269 6.419 2.58e-08 ***
## treatmentNP 0.01819 0.13269 0.137 0.89142
## treatmentP -0.26763 0.13916 -1.923 0.05929 .
## temperature12:treatmentN -0.26644 0.18765 -1.420 0.16089
## temperature16:treatmentN -0.64457 0.18765 -3.435 0.00109 **
## temperature12:treatmentNP -0.06016 0.18765 -0.321 0.74965
## temperature16:treatmentNP 0.16898 0.18765 0.901 0.37151
## temperature12:treatmentP 0.15352 0.19228 0.798 0.42782
## temperature16:treatmentP 0.19278 0.19228 1.003 0.32015
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2298 on 59 degrees of freedom
## Multiple R-squared: 0.9742, Adjusted R-squared: 0.9689
## F-statistic: 185.5 on 12 and 59 DF, p-value: < 2.2e-16
## Anova Table (Type III tests)
##
## Response: logNH4
## Sum Sq Df F value Pr(>F)
## temperature 33.739 3 212.9332 < 2.2e-16 ***
## treatment 4.058 3 25.6098 9.801e-11 ***
## temperature:treatment 1.466 6 4.6264 0.0006485 ***
## Residuals 3.116 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Call:
## lm(formula = logNH4 ~ temperature * treatment - 1, data = expchem)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.46007 -0.13249 -0.00846 0.12301 1.04024
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## temperature8 -1.39613 0.09382 -14.880 < 2e-16 ***
## temperature12 -1.23594 0.09382 -13.173 < 2e-16 ***
## temperature16 -1.46508 0.09382 -15.615 < 2e-16 ***
## treatmentN 0.85174 0.13269 6.419 2.58e-08 ***
## treatmentNP 0.01819 0.13269 0.137 0.89142
## treatmentP -0.26763 0.13916 -1.923 0.05929 .
## temperature12:treatmentN -0.26644 0.18765 -1.420 0.16089
## temperature16:treatmentN -0.64457 0.18765 -3.435 0.00109 **
## temperature12:treatmentNP -0.06016 0.18765 -0.321 0.74965
## temperature16:treatmentNP 0.16898 0.18765 0.901 0.37151
## temperature12:treatmentP 0.15352 0.19228 0.798 0.42782
## temperature16:treatmentP 0.19278 0.19228 1.003 0.32015
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2298 on 59 degrees of freedom
## Multiple R-squared: 0.9742, Adjusted R-squared: 0.9689
## F-statistic: 185.5 on 12 and 59 DF, p-value: < 2.2e-16
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 0.6553 0.7738
## 59
Pairwise comparsions for each level of temperature
## $emmeans
## temperature = 8:
## treatment emmean SE df lower.CL upper.CL
## C 0.04500000 0.02458718 59 -0.004198824 0.09419882
## N 0.28666667 0.02458718 59 0.237467842 0.33586549
## NP 0.04333333 0.02458718 59 -0.005865491 0.09253216
## P 0.02400000 0.02693390 59 -0.029894612 0.07789461
##
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C 0.07000000 0.02458718 59 0.020801176 0.11919882
## N 0.22666667 0.02458718 59 0.177467842 0.27586549
## NP 0.05666667 0.02458718 59 0.007467842 0.10586549
## P 0.10500000 0.02458718 59 0.055801176 0.15419882
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C 0.03500000 0.02458718 59 -0.014198824 0.08419882
## N 0.05833333 0.02458718 59 0.009134509 0.10753216
## NP 0.05666667 0.02458718 59 0.007467842 0.10586549
## P 0.03166667 0.02458718 59 -0.017532158 0.08086549
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 8:
## contrast estimate SE df t.ratio p.value
## C - N -0.241666667 0.03477152 59 -6.950 <.0001
## C - NP 0.001666667 0.03477152 59 0.048 1.0000
## C - P 0.021000000 0.03646867 59 0.576 0.9389
## N - NP 0.243333333 0.03477152 59 6.998 <.0001
## N - P 0.262666667 0.03646867 59 7.203 <.0001
## NP - P 0.019333333 0.03646867 59 0.530 0.9514
##
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N -0.156666667 0.03477152 59 -4.506 0.0002
## C - NP 0.013333333 0.03477152 59 0.383 0.9806
## C - P -0.035000000 0.03477152 59 -1.007 0.7462
## N - NP 0.170000000 0.03477152 59 4.889 <.0001
## N - P 0.121666667 0.03477152 59 3.499 0.0048
## NP - P -0.048333333 0.03477152 59 -1.390 0.5105
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N -0.023333333 0.03477152 59 -0.671 0.9076
## C - NP -0.021666667 0.03477152 59 -0.623 0.9243
## C - P 0.003333333 0.03477152 59 0.096 0.9997
## N - NP 0.001666667 0.03477152 59 0.048 1.0000
## N - P 0.026666667 0.03477152 59 0.767 0.8690
## NP - P 0.025000000 0.03477152 59 0.719 0.8892
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Pairwise comparsions for each level of nutrient treatment
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 0.04500000 0.02458718 59 -0.004198824 0.09419882
## 12 0.07000000 0.02458718 59 0.020801176 0.11919882
## 16 0.03500000 0.02458718 59 -0.014198824 0.08419882
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 0.28666667 0.02458718 59 0.237467842 0.33586549
## 12 0.22666667 0.02458718 59 0.177467842 0.27586549
## 16 0.05833333 0.02458718 59 0.009134509 0.10753216
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 0.04333333 0.02458718 59 -0.005865491 0.09253216
## 12 0.05666667 0.02458718 59 0.007467842 0.10586549
## 16 0.05666667 0.02458718 59 0.007467842 0.10586549
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 0.02400000 0.02693390 59 -0.029894612 0.07789461
## 12 0.10500000 0.02458718 59 0.055801176 0.15419882
## 16 0.03166667 0.02458718 59 -0.017532158 0.08086549
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -2.500000e-02 0.03477152 59 -0.719 0.7533
## 8 - 16 1.000000e-02 0.03477152 59 0.288 0.9555
## 12 - 16 3.500000e-02 0.03477152 59 1.007 0.5757
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 6.000000e-02 0.03477152 59 1.726 0.2043
## 8 - 16 2.283333e-01 0.03477152 59 6.567 <.0001
## 12 - 16 1.683333e-01 0.03477152 59 4.841 <.0001
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -1.333333e-02 0.03477152 59 -0.383 0.9223
## 8 - 16 -1.333333e-02 0.03477152 59 -0.383 0.9223
## 12 - 16 5.030698e-17 0.03477152 59 0.000 1.0000
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -8.100000e-02 0.03646867 59 -2.221 0.0759
## 8 - 16 -7.666667e-03 0.03646867 59 -0.210 0.9759
## 12 - 16 7.333333e-02 0.03477152 59 2.109 0.0967
##
## P value adjustment: tukey method for comparing a family of 3 estimates
## Error in UseMethod("cld"): no applicable method for 'cld' applied to an object of class "c('emm_list', 'list')"
## Error in UseMethod("cld"): no applicable method for 'cld' applied to an object of class "c('emm_list', 'list')"
Before we jump into any statistical analyses, just want to get an idea of how things play out by treatment and temperature .
| tempxtrt | treatment | temperature | meanchangeDOC_mgL | sdchangeDOC_mgL | n |
|---|---|---|---|---|---|
| 8Control | C | 8 | 2.759 | 0.627 | 6 |
| 8N+P | NP | 8 | 2.703 | 0.482 | 6 |
| 8Nitrogen | N | 8 | 3.136 | 1.338 | 6 |
| 8Phosphorus | P | 8 | 2.329 | 0.565 | 5 |
| 12Control | C | 12 | 1.898 | 0.347 | 6 |
| 12N+P | NP | 12 | 1.224 | 0.258 | 6 |
| 12Nitrogen | N | 12 | 1.771 | 0.226 | 6 |
| 12Phosphorus | P | 12 | 1.268 | 0.159 | 6 |
| 16Control | C | 16 | 2.534 | 0.344 | 6 |
| 16N+P | NP | 16 | 1.446 | 0.387 | 6 |
| 16Nitrogen | N | 16 | 2.841 | 0.481 | 6 |
| 16Phosphorus | P | 16 | 1.253 | 0.209 | 6 |
Here we have some summary statistics for DOC exudation. First we are looking at change in DOC from start to end of experiment by TEMPERTURE and by TREATMENT separately. I’ve also plotted boxplots with individual points. Due to the way I calculated these values (Initial DOC minus final DOC), negative numbers here mean there is MORE dissolved organic carbon at the end of the experiment than the begining.
## -------------------------------------------------------------------------
## changeDOC_mgL ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean 2.7e+00 1.5e+00 2.0e+00
## median 2.8e+00 1.4e+00 2.0e+00
## sd 8.3e-01 3.9e-01 7.8e-01
## IQR 9.5e-01 4.6e-01 1.3e+00
## n 2e+01 2e+01 2e+01
## np 32.4% 33.8% 33.8%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 24.495, df = 2, p-value = 4.796e-06
## -------------------------------------------------------------------------
## changeDOC_mgL ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N P NP
## mean 2.4e+00 2.6e+00 1.6e+00 1.8e+00
## median 2.3e+00 2.2e+00 1.4e+00 1.6e+00
## sd 5.7e-01 9.9e-01 5.9e-01 7.6e-01
## IQR 7.8e-01 1.2e+00 3.9e-01 1.1e+00
## n 2e+01 2e+01 2e+01 2e+01
## np 25.4% 25.4% 23.9% 25.4%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 20.608, df = 3, p-value = 0.000127
## # A tibble: 18 x 9
## temperature .y. group1 group2 p p.adj p.format p.signif method
## <fct> <chr> <chr> <chr> <dbl> <dbl> <chr> <chr> <chr>
## 1 8 chang… C N 0.936 1 0.9361 ns Wilco…
## 2 8 chang… C P 0.315 1 0.3153 ns Wilco…
## 3 8 chang… C NP 0.936 1 0.9362 ns Wilco…
## 4 8 chang… N P 0.234 1 0.2343 ns Wilco…
## 5 8 chang… N NP 0.810 1 0.8099 ns Wilco…
## 6 8 chang… P NP 0.583 1 0.5830 ns Wilco…
## 7 12 chang… C N 0.689 1 0.6889 ns Wilco…
## 8 12 chang… C P 0.00507 0.09 0.0051 ** Wilco…
## 9 12 chang… C NP 0.00630 0.09 0.0063 ** Wilco…
## 10 12 chang… N P 0.00507 0.09 0.0051 ** Wilco…
## 11 12 chang… N NP 0.00500 0.09 0.0050 ** Wilco…
## 12 12 chang… P NP 1 1 1.0000 ns Wilco…
## 13 16 chang… C N 0.261 1 0.2615 ns Wilco…
## 14 16 chang… C P 0.00507 0.09 0.0051 ** Wilco…
## 15 16 chang… C NP 0.00507 0.09 0.0051 ** Wilco…
## 16 16 chang… N P 0.00507 0.09 0.0051 ** Wilco…
## 17 16 chang… N NP 0.00507 0.09 0.0051 ** Wilco…
## 18 16 chang… P NP 0.423 1 0.4225 ns Wilco…
## $x
## [1] "Incubation temperature"
##
## $y
## [1] "Net change in DOC (mg/L)"
##
## attr(,"class")
## [1] "labels"
From our preliminary glance at the data, it looks like at lower temperatures, there is more DOC released (still pondering the ecological significance of that… fueling microbial loop?). At higher temperatures, it looks like there is less DOC at the end of the experiment.
## Anova Table (Type III tests)
##
## Response: changeDOC_mgL
## Sum Sq Df F value Pr(>F)
## (Intercept) 45.684 1 154.9900 < 2e-16 ***
## temperature 2.397 2 4.0659 0.02217 *
## treatment 1.794 3 2.0290 0.11959
## temperature:treatment 3.517 6 1.9887 0.08164 .
## Residuals 17.390 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), we see there is a evidence of a marginally significant interaction effect (F(6,59)=1.9887, p=0.08164). The test for the main effect of temperature treatment is significant (p=0.02217) but the test for the main effect of nutrient treatment alone is not significant (p=0.11959). After fitting an ANOVA model it is important to always check the relevant model assumptions. This includes making QQ-plots and residual plots. Will use Levene’s test of equality of variances too.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 2.0599 0.03806 *
## 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
It looks like there is just one potential outlier in the data. Overall, neither plot indicates a significant violation of normality assumptions. We’re good here.
## Anova Table (Type III tests)
##
## Response: changeDOC_mgL
## Sum Sq Df F value Pr(>F)
## (Intercept) 45.684 1 154.9900 < 2e-16 ***
## temperature 2.397 2 4.0659 0.02217 *
## treatment 1.794 3 2.0290 0.11959
## temperature:treatment 3.517 6 1.9887 0.08164 .
## Residuals 17.390 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 2.0599 0.03806 *
## 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
If we take the one super high DOC value out, we meet the assumptions of normality. Not sure if this is really ethical though.
Pairwise comparsions for each level of temperature. At 12 degrees, we see significnant differences between C - NP, N - P, and N - NP. At 16 degress, there is a difference between C - P, C-NP, N-P. There are no differences at 8 degrees.
## $emmeans
## temperature = 8:
## treatment emmean SE df lower.CL upper.CL
## C 2.759333 0.2216422 59 2.3158284 3.202838
## N 3.136000 0.2216422 59 2.6924951 3.579505
## P 2.329000 0.2427968 59 1.8431647 2.814835
## NP 2.702667 0.2216422 59 2.2591617 3.146172
##
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C 1.897667 0.2216422 59 1.4541617 2.341172
## N 1.771000 0.2216422 59 1.3274951 2.214505
## P 1.267667 0.2216422 59 0.8241617 1.711172
## NP 1.224333 0.2216422 59 0.7808284 1.667838
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C 2.534333 0.2216422 59 2.0908284 2.977838
## N 2.841000 0.2216422 59 2.3974951 3.284505
## P 1.253333 0.2216422 59 0.8098284 1.696838
## NP 1.446000 0.2216422 59 1.0024951 1.889505
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 8:
## contrast estimate SE df t.ratio p.value
## C - N -0.37666667 0.3134493 59 -1.202 0.6282
## C - P 0.43033333 0.3287484 59 1.309 0.5610
## C - NP 0.05666667 0.3134493 59 0.181 0.9979
## N - P 0.80700000 0.3287484 59 2.455 0.0780
## N - NP 0.43333333 0.3134493 59 1.382 0.5152
## P - NP -0.37366667 0.3287484 59 -1.137 0.6686
##
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N 0.12666667 0.3134493 59 0.404 0.9775
## C - P 0.63000000 0.3134493 59 2.010 0.1960
## C - NP 0.67333333 0.3134493 59 2.148 0.1501
## N - P 0.50333333 0.3134493 59 1.606 0.3833
## N - NP 0.54666667 0.3134493 59 1.744 0.3106
## P - NP 0.04333333 0.3134493 59 0.138 0.9990
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N -0.30666667 0.3134493 59 -0.978 0.7622
## C - P 1.28100000 0.3134493 59 4.087 0.0008
## C - NP 1.08833333 0.3134493 59 3.472 0.0052
## N - P 1.58766667 0.3134493 59 5.065 <.0001
## N - NP 1.39500000 0.3134493 59 4.450 0.0002
## P - NP -0.19266667 0.3134493 59 -0.615 0.9270
##
## P value adjustment: tukey method for comparing a family of 4 estimates
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 2.759333 0.2216422 59 2.3158284 3.202838
## 12 1.897667 0.2216422 59 1.4541617 2.341172
## 16 2.534333 0.2216422 59 2.0908284 2.977838
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 3.136000 0.2216422 59 2.6924951 3.579505
## 12 1.771000 0.2216422 59 1.3274951 2.214505
## 16 2.841000 0.2216422 59 2.3974951 3.284505
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 2.329000 0.2427968 59 1.8431647 2.814835
## 12 1.267667 0.2216422 59 0.8241617 1.711172
## 16 1.253333 0.2216422 59 0.8098284 1.696838
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 2.702667 0.2216422 59 2.2591617 3.146172
## 12 1.224333 0.2216422 59 0.7808284 1.667838
## 16 1.446000 0.2216422 59 1.0024951 1.889505
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 0.86166667 0.3134493 59 2.749 0.0213
## 8 - 16 0.22500000 0.3134493 59 0.718 0.7540
## 12 - 16 -0.63666667 0.3134493 59 -2.031 0.1137
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 1.36500000 0.3134493 59 4.355 0.0002
## 8 - 16 0.29500000 0.3134493 59 0.941 0.6167
## 12 - 16 -1.07000000 0.3134493 59 -3.414 0.0033
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 1.06133333 0.3287484 59 3.228 0.0057
## 8 - 16 1.07566667 0.3287484 59 3.272 0.0050
## 12 - 16 0.01433333 0.3134493 59 0.046 0.9988
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 1.47833333 0.3134493 59 4.716 <.0001
## 8 - 16 1.25666667 0.3134493 59 4.009 0.0005
## 12 - 16 -0.22166667 0.3134493 59 -0.707 0.7602
##
## P value adjustment: tukey method for comparing a family of 3 estimates
Pairwise comparsions for each level of nutrient treatment. In the control, there is a difference between 12-16 degrees. In the P treatment, there is a difference between 8-12, and 8-16. And in the NP treatment, there is a difference between 8-12, and 8-16.
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 2.759333 0.2216422 59 2.3158284 3.202838
## 12 1.897667 0.2216422 59 1.4541617 2.341172
## 16 2.534333 0.2216422 59 2.0908284 2.977838
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 3.136000 0.2216422 59 2.6924951 3.579505
## 12 1.771000 0.2216422 59 1.3274951 2.214505
## 16 2.841000 0.2216422 59 2.3974951 3.284505
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 2.329000 0.2427968 59 1.8431647 2.814835
## 12 1.267667 0.2216422 59 0.8241617 1.711172
## 16 1.253333 0.2216422 59 0.8098284 1.696838
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 2.702667 0.2216422 59 2.2591617 3.146172
## 12 1.224333 0.2216422 59 0.7808284 1.667838
## 16 1.446000 0.2216422 59 1.0024951 1.889505
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 0.86166667 0.3134493 59 2.749 0.0213
## 8 - 16 0.22500000 0.3134493 59 0.718 0.7540
## 12 - 16 -0.63666667 0.3134493 59 -2.031 0.1137
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 1.36500000 0.3134493 59 4.355 0.0002
## 8 - 16 0.29500000 0.3134493 59 0.941 0.6167
## 12 - 16 -1.07000000 0.3134493 59 -3.414 0.0033
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 1.06133333 0.3287484 59 3.228 0.0057
## 8 - 16 1.07566667 0.3287484 59 3.272 0.0050
## 12 - 16 0.01433333 0.3134493 59 0.046 0.9988
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 1.47833333 0.3134493 59 4.716 <.0001
## 8 - 16 1.25666667 0.3134493 59 4.009 0.0005
## 12 - 16 -0.22166667 0.3134493 59 -0.707 0.7602
##
## P value adjustment: tukey method for comparing a family of 3 estimates
DOCrelease_summary <- expchem %>%
group_by(temperature) %>%
summarize(meanchangeDOC_mgL=mean(changeDOC_mgL),
sdchangeDOC_mgL=sd(changeDOC_mgL),
n= n()) %>%
arrange(temperature)
DOCrelease_summary
## # A tibble: 3 x 4
## temperature meanchangeDOC_mgL sdchangeDOC_mgL n
## <fct> <dbl> <dbl> <int>
## 1 8 2.75 0.830 23
## 2 12 1.54 0.387 24
## 3 16 2.02 0.776 24
compare_means(changeDOC_mgL ~ temperature , data = expchem_trim)
## # A tibble: 3 x 8
## .y. group1 group2 p p.adj p.format p.signif method
## <chr> <chr> <chr> <dbl> <dbl> <chr> <chr> <chr>
## 1 changeDOC_… 8 12 4.55e-7 0.0000014 4.5e-07 **** Wilcox…
## 2 changeDOC_… 8 16 6.24e-3 0.012 0.0062 ** Wilcox…
## 3 changeDOC_… 12 16 5.78e-2 0.058 0.0578 ns Wilcox…
temp_DOC_comparisons <- list(c("8","12"),c("8","16"))
DOCrelease_temp_boxplot <- ggplot(expchem_trim, aes(x=temperature, y=changeDOC_mgL, fill=temperature)) +
geom_boxplot()+
geom_point(size=3, alpha=0.4,position = position_jitter()) +
stat_boxplot(geom ='errorbar') +
# scale_fill_manual(values=trtcolors,
# name="Nutrient\nTreatment",
# labels=c("Control", "Nitrogen", "Phosphorus", "N&P"))+
coord_cartesian(ylim=c(0,5))+
scale_y_continuous(breaks=seq(0, 5, 1), name="Net change DOC (mg/L)")+
theme(axis.text.x = element_text(colour = "black", size = 24), axis.title.y=element_text(size=24, vjust=1.2, face="bold"),
axis.text.y = element_text(colour = "black", size = 24), axis.title.x=element_text(size=24, vjust=-0.1, face="bold"),
legend.text = element_text(colour="black", size = 24), legend.title=element_text(size=24, face="bold"),
legend.key.size = unit(3, 'lines'),
legend.position="none",
plot.title=element_text(size=24,face="bold"),
panel.border = element_rect(size=1.0, fill=NA),
panel.grid.minor = element_blank(),
panel.grid.major = element_blank())+
stat_compare_means(comparisons=temp_DOC_comparisons, label="p.signif")+
stat_compare_means(label.y=5)+
labs(x="Temperature")
DOCrelease_temp_boxplot
# ggsave("DOCchange_bytemp.png", width=9.5, height=5,units="in")
Before we jump into any statistical analyses, just want to get an idea of how things play out by treatment and temperature .
| tempxtrt | treatment | temperature | meanchangePO4_uM | sdchangePO4_uM | n |
|---|---|---|---|---|---|
| 8Control | C | 8 | -5.455 | 2.146 | 6 |
| 8N+P | NP | 8 | -1.083 | 1.103 | 6 |
| 8Nitrogen | N | 8 | -8.886 | 1.358 | 6 |
| 8Phosphorus | P | 8 | 1.668 | 0.267 | 5 |
| 12Control | C | 12 | -4.806 | 2.353 | 6 |
| 12N+P | NP | 12 | 7.893 | 0.000 | 6 |
| 12Nitrogen | N | 12 | -5.016 | 4.002 | 6 |
| 12Phosphorus | P | 12 | 7.493 | 0.000 | 6 |
| 16Control | C | 16 | -0.007 | 0.000 | 6 |
| 16N+P | NP | 16 | 7.893 | 0.000 | 6 |
| 16Nitrogen | N | 16 | -0.007 | 0.000 | 6 |
| 16Phosphorus | P | 16 | 7.493 | 0.000 | 6 |
Here we have some summary statistics for PO4 uptake/release. First we are looking at change in PO4 from start to end of experiment by TEMPERTURE and by TREATMENT separately. I’ve also plotted boxplots with individual points. Due to the way I calculated these values (Initial PO4 minus final PO4), negative numbers here mean there is MORE phosphate at the end of the experiment than the begining. Positive number indiciate uptake… slightly counterintuitive.
## -------------------------------------------------------------------------
## changePO4_uM ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean -3.7e+00 1.4e+00 3.8e+00
## median -3.3e+00 3.7e+00 3.7e+00
## sd 4.3e+00 6.8e+00 3.9e+00
## IQR 6.9e+00 1.3e+01 7.6e+00
## n 2e+01 2e+01 2e+01
## np 32.4% 33.8% 33.8%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 20.651, df = 2, p-value = 3.279e-05
## -------------------------------------------------------------------------
## changePO4_uM ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N P NP
## mean -3.4e+00 -4.6e+00 5.8e+00 4.9e+00
## median -5.7e+00 -6.7e+00 7.5e+00 7.9e+00
## sd 3.0e+00 4.4e+00 2.7e+00 4.4e+00
## IQR 5.9e+00 8.4e+00 5.5e+00 8.3e+00
## n 2e+01 2e+01 2e+01 2e+01
## np 25.4% 25.4% 23.9% 25.4%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 39.501, df = 3, p-value = 1.359e-08
From our preliminary glance at the data, it looks like PO4 was released (averaging over all nutrient treatments). Interestingly, in the P & NP treatments, there was PO4 uptake compared to C & N treatments, where there was more PO4 leftover in the begining of the experiment than there were at the beginning.
## Anova Table (Type III tests)
##
## Response: changePO4_uM
## Sum Sq Df F value Pr(>F)
## (Intercept) 178.54 1 71.9709 8.366e-12 ***
## temperature 106.27 2 21.4189 1.016e-07 ***
## treatment 364.96 3 49.0382 4.970e-16 ***
## temperature:treatment 143.09 6 9.6131 2.305e-07 ***
## Residuals 146.36 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), we see there is a evidence of a marginally significant interaction effect (F(6,59)=9.6131, p<0.0001). The tests for the main effects of temperature and nutrient treatments are BOTH significant. After fitting an ANOVA model it is important to always check the relevant model assumptions. This includes making QQ-plots and residual plots. Will use Levene’s test of equality of variances too.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 2.9936 0.003176 **
## 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
PO4 uptake does not pass the Levene’s test and the Normal Q-Q plot doesn’t look good. We can try log transforming the data and see if that helps.
## Anova Table (Type III tests)
##
## Response: log_changePO4_uM
## Sum Sq Df F value Pr(>F)
## (Intercept) 4.2420 1 353.0388 < 2.2e-16 ***
## temperature 0.2352 2 9.7859 0.0002137 ***
## treatment 1.3423 3 37.2365 1.239e-13 ***
## temperature:treatment 0.4092 6 5.6765 0.0001046 ***
## Residuals 0.7089 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 2.5357 0.01079 *
## 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
This is slightly better, but there is still quite a bit spread in our residuals which is making us fail the Levene test. Might need to think of a workaround here.
Pairwise comparsions for each level of temperature. Pretty much the only time we don’t see a significant different between treatments is at 8 degree between NP & P, 12 degrees between C & N and P & NP, and 16 degrees between C & N, and P & NP.
## $emmeans
## temperature = 8:
## treatment emmean SE df lower.CL upper.CL
## C 0.8408372 0.04475079 59 0.7512911 0.9303833
## N 0.5230641 0.04475079 59 0.4335180 0.6126103
## P 1.1541683 0.04902204 59 1.0560754 1.2522611
## NP 1.0593662 0.04475079 59 0.9698201 1.1489123
##
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C 0.8786424 0.04475079 59 0.7890962 0.9681885
## N 0.8308817 0.04475079 59 0.7413356 0.9204278
## P 1.3029583 0.04475079 59 1.2134122 1.3925044
## NP 1.3115208 0.04475079 59 1.2219746 1.4010669
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C 1.0999912 0.04475079 59 1.0104451 1.1895374
## N 1.0999912 0.04475079 59 1.0104451 1.1895374
## P 1.3029583 0.04475079 59 1.2134122 1.3925044
## NP 1.3115208 0.04475079 59 1.2219746 1.4010669
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 8:
## contrast estimate SE df t.ratio p.value
## C - N 3.177731e-01 0.06328718 59 5.021 <.0001
## C - P -3.133311e-01 0.06637615 59 -4.721 0.0001
## C - NP -2.185290e-01 0.06328718 59 -3.453 0.0055
## N - P -6.311041e-01 0.06637615 59 -9.508 <.0001
## N - NP -5.363021e-01 0.06328718 59 -8.474 <.0001
## P - NP 9.480208e-02 0.06637615 59 1.428 0.4871
##
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N 4.776069e-02 0.06328718 59 0.755 0.8743
## C - P -4.243159e-01 0.06328718 59 -6.705 <.0001
## C - NP -4.328784e-01 0.06328718 59 -6.840 <.0001
## N - P -4.720766e-01 0.06328718 59 -7.459 <.0001
## N - NP -4.806391e-01 0.06328718 59 -7.595 <.0001
## P - NP -8.562444e-03 0.06328718 59 -0.135 0.9991
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N -2.775558e-16 0.06328718 59 0.000 1.0000
## C - P -2.029671e-01 0.06328718 59 -3.207 0.0113
## C - NP -2.115295e-01 0.06328718 59 -3.342 0.0077
## N - P -2.029671e-01 0.06328718 59 -3.207 0.0113
## N - NP -2.115295e-01 0.06328718 59 -3.342 0.0077
## P - NP -8.562444e-03 0.06328718 59 -0.135 0.9991
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Pairwise comparsions for each level of nutrient treatment. In the control nutrient treatment, there is a different between 8 and 16 degrees and 12 and 16 degrees. In the nitrogen nutrient treatment, there is a statistically significant difference with temperature. In the P treatment, there is a statastically significant difference between 8 and 12 and 8 and 16, but not between 12 and 16. Finally, in the NP treatment, therre is a difference between 8 and 12 and 8 and 16, but not between 12 and 16.
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 0.8408372 0.04475079 59 0.7512911 0.9303833
## 12 0.8786424 0.04475079 59 0.7890962 0.9681885
## 16 1.0999912 0.04475079 59 1.0104451 1.1895374
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 0.5230641 0.04475079 59 0.4335180 0.6126103
## 12 0.8308817 0.04475079 59 0.7413356 0.9204278
## 16 1.0999912 0.04475079 59 1.0104451 1.1895374
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 1.1541683 0.04902204 59 1.0560754 1.2522611
## 12 1.3029583 0.04475079 59 1.2134122 1.3925044
## 16 1.3029583 0.04475079 59 1.2134122 1.3925044
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 1.0593662 0.04475079 59 0.9698201 1.1489123
## 12 1.3115208 0.04475079 59 1.2219746 1.4010669
## 16 1.3115208 0.04475079 59 1.2219746 1.4010669
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -3.780519e-02 0.06328718 59 -0.597 0.8221
## 8 - 16 -2.591540e-01 0.06328718 59 -4.095 0.0004
## 12 - 16 -2.213489e-01 0.06328718 59 -3.498 0.0026
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -3.078176e-01 0.06328718 59 -4.864 <.0001
## 8 - 16 -5.769271e-01 0.06328718 59 -9.116 <.0001
## 12 - 16 -2.691095e-01 0.06328718 59 -4.252 0.0002
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -1.487901e-01 0.06637615 59 -2.242 0.0725
## 8 - 16 -1.487901e-01 0.06637615 59 -2.242 0.0725
## 12 - 16 -9.714451e-17 0.06328718 59 0.000 1.0000
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -2.521546e-01 0.06328718 59 -3.984 0.0005
## 8 - 16 -2.521546e-01 0.06328718 59 -3.984 0.0005
## 12 - 16 -3.053113e-16 0.06328718 59 0.000 1.0000
##
## P value adjustment: tukey method for comparing a family of 3 estimates
Tutorial found here: http://rcompanion.org/rcompanion/d_08a.html More info on non-parametric tests and R packages here: https://journal.r-project.org/archive/2016/RJ-2016-027/RJ-2016-027.pdf
First, Produce Huber M-estimators and confidence intervals by group NOT PRINTING AT THIS TIME
Here we have some summary statistics for PO4 uptake/release. First we are looking at change in PO4 from start to end of experiment by TEMPERTURE and by TREATMENT separately. I’ve also plotted boxplots with individual points. Due to the way I calculated these values (Initial PO4 minus final PO4), negative numbers here mean there is MORE phosphate at the end of the experiment than the begining. Positive number indiciate uptake… slightly counterintuitive.
## -------------------------------------------------------------------------
## PO4_uM ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean 8.2e+00 3.3e+00 8.1e-01
## median 8.6e+00 8.1e-01 8.1e-01
## sd 2.2e+00 3.3e+00 0.0
## IQR 2.5e+00 5.7e+00 0.0
## n 2e+01 2e+01 2e+01
## np 32.4% 33.8% 33.8%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 48.704, df = 2, p-value = 2.656e-11
## -------------------------------------------------------------------------
## PO4_uM ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N NP P
## mean 4.2e+00 5.4e+00 3.8e+00 2.5e+00
## median 6.5e+00 7.5e+00 8.1e-01 8.1e-01
## sd 3.0e+00 4.4e+00 4.4e+00 2.7e+00
## IQR 5.9e+00 8.4e+00 8.3e+00 5.5e+00
## n 2e+01 2e+01 2e+01 2e+01
## np 25.4% 25.4% 25.4% 23.9%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 5.4873, df = 3, p-value = 0.1394
Before we jump into any statistical analyses, just want to get an idea of how things play out by treatment and temperature .
| tempxtrt | treatment | temperature | meanphytochla | sdphytochla | n |
|---|---|---|---|---|---|
| 8Control | C | 8 | 19.947 | 6.011 | 6 |
| 8N+P | NP | 8 | 71.782 | 5.824 | 6 |
| 8Nitrogen | N | 8 | 19.908 | 6.649 | 6 |
| 8Phosphorus | P | 8 | 27.330 | 5.930 | 5 |
| 12Control | C | 12 | 22.712 | 7.839 | 6 |
| 12N+P | NP | 12 | 65.117 | 5.632 | 6 |
| 12Nitrogen | N | 12 | 21.140 | 8.846 | 6 |
| 12Phosphorus | P | 12 | 17.715 | 6.384 | 6 |
| 16Control | C | 16 | 19.338 | 7.403 | 6 |
| 16N+P | NP | 16 | 42.807 | 12.355 | 6 |
| 16Nitrogen | N | 16 | 37.978 | 6.251 | 6 |
| 16Phosphorus | P | 16 | 19.847 | 4.565 | 6 |
Here we have some summary statistics for Beaker Chla. First we are looking at how much chlorophyll a there was in the beaker at the end of experiment by TEMPERTURE and by TREATMENT separately. I’ve also plotted boxplots with individual points. Not sure where all these phytoplankton came from, since we used filtered lake water.
## -------------------------------------------------------------------------
## chl_a_phyto ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean 3.5e+01 3.2e+01 3.0e+01
## median 2.5e+01 2.3e+01 2.9e+01
## sd 2.3e+01 2.1e+01 1.3e+01
## IQR 3.2e+01 2.4e+01 1.9e+01
## n 2e+01 2e+01 2e+01
## np 32.4% 33.8% 33.8%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 0.37119, df = 2, p-value = 0.8306
## -------------------------------------------------------------------------
## chl_a_phyto ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N P NP
## mean 2.1e+01 2.6e+01 2.1e+01 6.0e+01
## median 2.2e+01 2.8e+01 2.2e+01 6.4e+01
## sd 6.9e+00 1.1e+01 6.7e+00 1.5e+01
## IQR 8.0e+00 1.9e+01 8.2e+00 1.6e+01
## n 2e+01 2e+01 2e+01 2e+01
## np 25.4% 25.4% 23.9% 25.4%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 35.38, df = 3, p-value = 1.013e-07
From our preliminary glance at the data, it doesn’t look like there is much a difference in the means, although in the 8 degree treatments, the spread was quite wide. Comparing across nutrient treatments, NP has the highest mean chla with N as the runner up.
## Anova Table (Type III tests)
##
## Response: chl_a_phyto
## Sum Sq Df F value Pr(>F)
## (Intercept) 2387.2 1 45.292 7.740e-09 ***
## temperature 38.8 2 0.368 0.6937
## treatment 11137.6 3 70.436 < 2.2e-16 ***
## temperature:treatment 4042.5 6 12.783 3.445e-09 ***
## Residuals 3109.8 59
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), we see there is a evidence of a significant interaction effect (F(6,59)=12.783, p<0.0001). The tests for the main effects of temperature is not significant, but the effect ot treatment is highly significant (F(3,59)=70.436, p<0.0001).
After fitting an ANOVA model it is important to always check the relevant model assumptions. This includes making QQ-plots and residual plots. Will use Levene’s test of equality of variances too.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 0.3043 0.9824
## 59
Data fit the model assumptions. We can move on.
Pairwise comparsions for each level of temperature. At 8 degrees, NP is different than C, N, and P. At 12 degrees, NP is different than C, N, and P. At 16 degrees, both N and NP are significantly different than the Control. N vs P and P versus NP are both significantly different from each other, but N and NP are not significantly different.
## $emmeans
## temperature = 8:
## treatment emmean SE df lower.CL upper.CL
## C 19.94667 2.963885 59 14.01595 25.87739
## N 19.90833 2.963885 59 13.97761 25.83905
## P 27.33000 3.246773 59 20.83322 33.82678
## NP 71.78167 2.963885 59 65.85095 77.71239
##
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C 22.71167 2.963885 59 16.78095 28.64239
## N 21.14000 2.963885 59 15.20928 27.07072
## P 17.71500 2.963885 59 11.78428 23.64572
## NP 65.11667 2.963885 59 59.18595 71.04739
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C 19.33833 2.963885 59 13.40761 25.26905
## N 37.97833 2.963885 59 32.04761 43.90905
## P 19.84667 2.963885 59 13.91595 25.77739
## NP 42.80667 2.963885 59 36.87595 48.73739
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 8:
## contrast estimate SE df t.ratio p.value
## C - N 0.03833333 4.191566 59 0.009 1.0000
## C - P -7.38333333 4.396151 59 -1.679 0.3435
## C - NP -51.83500000 4.191566 59 -12.367 <.0001
## N - P -7.42166667 4.396151 59 -1.688 0.3390
## N - NP -51.87333333 4.191566 59 -12.376 <.0001
## P - NP -44.45166667 4.396151 59 -10.111 <.0001
##
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N 1.57166667 4.191566 59 0.375 0.9819
## C - P 4.99666667 4.191566 59 1.192 0.6342
## C - NP -42.40500000 4.191566 59 -10.117 <.0001
## N - P 3.42500000 4.191566 59 0.817 0.8461
## N - NP -43.97666667 4.191566 59 -10.492 <.0001
## P - NP -47.40166667 4.191566 59 -11.309 <.0001
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N -18.64000000 4.191566 59 -4.447 0.0002
## C - P -0.50833333 4.191566 59 -0.121 0.9994
## C - NP -23.46833333 4.191566 59 -5.599 <.0001
## N - P 18.13166667 4.191566 59 4.326 0.0003
## N - NP -4.82833333 4.191566 59 -1.152 0.6591
## P - NP -22.96000000 4.191566 59 -5.478 <.0001
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Pairwise comparsions for each level of nutrient treatment. In the control nutrient treatment, there is a different between 8 and 16 degrees and 12 and 16 degrees. In the nitrogen nutrient treatment, there is a statistically significant difference with temperature. In the P treatment, there is a statastically significant difference between 8 and 12 and 8 and 16, but not between 12 and 16. Finally, in the NP treatment, therre is a difference between 8 and 12 and 8 and 16, but not between 12 and 16.
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 0.8408372 0.04475079 59 0.7512911 0.9303833
## 12 0.8786424 0.04475079 59 0.7890962 0.9681885
## 16 1.0999912 0.04475079 59 1.0104451 1.1895374
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 0.5230641 0.04475079 59 0.4335180 0.6126103
## 12 0.8308817 0.04475079 59 0.7413356 0.9204278
## 16 1.0999912 0.04475079 59 1.0104451 1.1895374
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 1.1541683 0.04902204 59 1.0560754 1.2522611
## 12 1.3029583 0.04475079 59 1.2134122 1.3925044
## 16 1.3029583 0.04475079 59 1.2134122 1.3925044
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 1.0593662 0.04475079 59 0.9698201 1.1489123
## 12 1.3115208 0.04475079 59 1.2219746 1.4010669
## 16 1.3115208 0.04475079 59 1.2219746 1.4010669
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -3.780519e-02 0.06328718 59 -0.597 0.8221
## 8 - 16 -2.591540e-01 0.06328718 59 -4.095 0.0004
## 12 - 16 -2.213489e-01 0.06328718 59 -3.498 0.0026
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -3.078176e-01 0.06328718 59 -4.864 <.0001
## 8 - 16 -5.769271e-01 0.06328718 59 -9.116 <.0001
## 12 - 16 -2.691095e-01 0.06328718 59 -4.252 0.0002
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -1.487901e-01 0.06637615 59 -2.242 0.0725
## 8 - 16 -1.487901e-01 0.06637615 59 -2.242 0.0725
## 12 - 16 -9.714451e-17 0.06328718 59 0.000 1.0000
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -2.521546e-01 0.06328718 59 -3.984 0.0005
## 8 - 16 -2.521546e-01 0.06328718 59 -3.984 0.0005
## 12 - 16 -3.053113e-16 0.06328718 59 0.000 1.0000
##
## P value adjustment: tukey method for comparing a family of 3 estimates
Before we jump into any statistical analyses, just want to get an idea of how things play out by treatment and temperature .
| treatment | temperature | meanGPP | sdGPP | n |
|---|---|---|---|---|
| C | 8 | 41.000 | 16.646 | 6 |
| N | 8 | 38.483 | 6.845 | 6 |
| NP | 8 | 45.567 | 19.708 | 6 |
| P | 8 | 51.250 | 6.149 | 6 |
| C | 12 | 60.517 | 17.748 | 6 |
| N | 12 | 38.167 | 6.124 | 6 |
| NP | 12 | 37.383 | 8.541 | 6 |
| P | 12 | 40.467 | 5.545 | 6 |
| C | 16 | 61.100 | 4.330 | 6 |
| N | 16 | 54.300 | 4.951 | 6 |
| NP | 16 | 58.217 | 6.065 | 6 |
| P | 16 | 52.483 | 3.459 | 6 |
Here we have some summary statistics for GPP. First we are looking at GPP by TEMPERTURE and by TREATMENT separately.
## -------------------------------------------------------------------------
## GPP_ugcm2h ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean 4.4e+01 4.4e+01 5.7e+01
## median 4.3e+01 4.3e+01 5.8e+01
## sd 1.4e+01 1.4e+01 5.6e+00
## IQR 2.0e+01 1.1e+01 9.6e+00
## n 2e+01 2e+01 2e+01
## np 33.3% 33.3% 33.3%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 22.444, df = 2, p-value = 1.337e-05
## -------------------------------------------------------------------------
## GPP_ugcm2h ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N P NP
## mean 5.4e+01 4.4e+01 4.8e+01 4.7e+01
## median 5.7e+01 4.3e+01 4.9e+01 4.7e+01
## sd 1.6e+01 9.6e+00 7.4e+00 1.5e+01
## IQR 1.9e+01 1.2e+01 8.0e+00 2.8e+01
## n 2e+01 2e+01 2e+01 2e+01
## np 25.0% 25.0% 25.0% 25.0%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 5.5192, df = 3, p-value = 0.1375
From our preliminary glance at the data, it looks like GPP might be higher on average in the warmest temperature treatment (16 degrees). As for the nutrients, it the controls have the widest spread in GPP and as a result have the highest mean.
## Anova Table (Type III tests)
##
## Response: GPP_ugcm2h
## Sum Sq Df F value Pr(>F)
## temperature 54459 3 167.7790 < 2e-16 ***
## treatment 567 3 1.7455 0.16733
## temperature:treatment 1970 6 3.0340 0.01171 *
## Residuals 6492 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), we see there is a evidence of a significant interaction effect (F(6,60)=3.0340, p=0.011710). The tests for the main effects of treatment is not significant, but the effect of temperature is significant (F(2,60)=7.2577, p=0.001504).
After fitting an ANOVA model it is important to always check the relevant model assumptions. This includes making QQ-plots and residual plots. Will use Levene’s test of equality of variances too.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 2.5787 0.009477 **
## 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Data doesn’t fit model assumptions, unfortunately. Let’s try a log transformation…
## Anova Table (Type III tests)
##
## Response: log_GPP_ugcm2h
## Sum Sq Df F value Pr(>F)
## (Intercept) 7.4772 1 108.8152 4.215e-15 ***
## temperature 0.8033 2 5.8450 0.004795 **
## treatment 0.4700 3 2.2800 0.088469 .
## temperature:treatment 0.8508 6 2.0636 0.070960 .
## Residuals 4.1229 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 3.6164 0.0005946 ***
## 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Anova Table (Type III tests)
##
## Response: sqrt_GPP
## Sum Sq Df F value Pr(>F)
## (Intercept) 238.241 1 423.1044 < 2.2e-16 ***
## temperature 8.570 2 7.6095 0.001134 **
## treatment 3.359 3 1.9887 0.125263
## temperature:treatment 10.324 6 3.0557 0.011245 *
## Residuals 33.785 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 2.9748 0.003271 **
## 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Neither log or square root transformations help. Will need to try something different. In the mean time, let’s look at some pairwise comparisons ANYWAY!
Pairwise comparsions for each level of temperature. At 8 degrees, there is no difference between treatments. At 12 degrees, the Control is different than all the other nutrient treatments. And at 16 degrees there is no different between treatments.
## $emmeans
## temperature = 8:
## treatment emmean SE df lower.CL upper.CL
## C 41.00000 4.246481 60 32.50577 49.49423
## N 38.48333 4.246481 60 29.98911 46.97756
## P 51.25000 4.246481 60 42.75577 59.74423
## NP 45.56667 4.246481 60 37.07244 54.06089
##
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C 60.51667 4.246481 60 52.02244 69.01089
## N 38.16667 4.246481 60 29.67244 46.66089
## P 40.46667 4.246481 60 31.97244 48.96089
## NP 37.38333 4.246481 60 28.88911 45.87756
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C 61.10000 4.246481 60 52.60577 69.59423
## N 54.30000 4.246481 60 45.80577 62.79423
## P 52.48333 4.246481 60 43.98911 60.97756
## NP 58.21667 4.246481 60 49.72244 66.71089
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 8:
## contrast estimate SE df t.ratio p.value
## C - N 2.5166667 6.005431 60 0.419 0.9750
## C - P -10.2500000 6.005431 60 -1.707 0.3292
## C - NP -4.5666667 6.005431 60 -0.760 0.8718
## N - P -12.7666667 6.005431 60 -2.126 0.1567
## N - NP -7.0833333 6.005431 60 -1.179 0.6420
## P - NP 5.6833333 6.005431 60 0.946 0.7800
##
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N 22.3500000 6.005431 60 3.722 0.0024
## C - P 20.0500000 6.005431 60 3.339 0.0077
## C - NP 23.1333333 6.005431 60 3.852 0.0016
## N - P -2.3000000 6.005431 60 -0.383 0.9807
## N - NP 0.7833333 6.005431 60 0.130 0.9992
## P - NP 3.0833333 6.005431 60 0.513 0.9556
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N 6.8000000 6.005431 60 1.132 0.6712
## C - P 8.6166667 6.005431 60 1.435 0.4830
## C - NP 2.8833333 6.005431 60 0.480 0.9632
## N - P 1.8166667 6.005431 60 0.303 0.9903
## N - NP -3.9166667 6.005431 60 -0.652 0.9144
## P - NP -5.7333333 6.005431 60 -0.955 0.7754
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Pairwise comparsions for each level of nutrient treatment. In the control treatment, there is a difference between 8 and 12 and 8 and 16, but none between 12 and 16. In the nitrogen treatment, there is a different between 8 and 16 ad 12 and 16, but not between 8 and 12. In the P treatment, there are no differences between temperatures. In the NP treatment, there is a difference between 8 and 16 and 12 and 16 but none between 8 and 12.
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 41.00000 4.246481 60 32.50577 49.49423
## 12 60.51667 4.246481 60 52.02244 69.01089
## 16 61.10000 4.246481 60 52.60577 69.59423
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 38.48333 4.246481 60 29.98911 46.97756
## 12 38.16667 4.246481 60 29.67244 46.66089
## 16 54.30000 4.246481 60 45.80577 62.79423
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 51.25000 4.246481 60 42.75577 59.74423
## 12 40.46667 4.246481 60 31.97244 48.96089
## 16 52.48333 4.246481 60 43.98911 60.97756
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 45.56667 4.246481 60 37.07244 54.06089
## 12 37.38333 4.246481 60 28.88911 45.87756
## 16 58.21667 4.246481 60 49.72244 66.71089
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -19.5166667 6.005431 60 -3.250 0.0053
## 8 - 16 -20.1000000 6.005431 60 -3.347 0.0040
## 12 - 16 -0.5833333 6.005431 60 -0.097 0.9948
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 0.3166667 6.005431 60 0.053 0.9985
## 8 - 16 -15.8166667 6.005431 60 -2.634 0.0285
## 12 - 16 -16.1333333 6.005431 60 -2.686 0.0249
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 10.7833333 6.005431 60 1.796 0.1797
## 8 - 16 -1.2333333 6.005431 60 -0.205 0.9770
## 12 - 16 -12.0166667 6.005431 60 -2.001 0.1207
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 8.1833333 6.005431 60 1.363 0.3668
## 8 - 16 -12.6500000 6.005431 60 -2.106 0.0971
## 12 - 16 -20.8333333 6.005431 60 -3.469 0.0028
##
## P value adjustment: tukey method for comparing a family of 3 estimates
Before we jump into any statistical analyses, just want to get an idea of how things play out by treatment and temperature .
| treatment | temperature | meanER | sdER | n |
|---|---|---|---|---|
| C | 8 | 0.517 | 4.449 | 6 |
| N | 8 | -13.533 | 5.007 | 6 |
| NP | 8 | -21.333 | 10.945 | 6 |
| P | 8 | -13.983 | 3.071 | 6 |
| C | 12 | -26.450 | 12.501 | 6 |
| N | 12 | -14.333 | 4.717 | 6 |
| NP | 12 | -10.783 | 5.969 | 6 |
| P | 12 | -18.783 | 6.492 | 6 |
| C | 16 | -42.067 | 2.196 | 6 |
| N | 16 | -35.217 | 4.990 | 6 |
| NP | 16 | -39.783 | 6.218 | 6 |
| P | 16 | -42.583 | 5.312 | 6 |
Here we have some summary statistics for ER. First we are looking at ER by TEMPERTURE and by TREATMENT separately.
## -------------------------------------------------------------------------
## ER_ugcm2h ~ temperature
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 3
##
##
## 8 12 16
## mean -1.2e+01 -1.8e+01 -4.0e+01
## median -1.3e+01 -1.5e+01 -4.1e+01
## sd 1.0e+01 9.6e+00 5.5e+00
## IQR 9.2e+00 9.6e+00 6.6e+00
## n 2e+01 2e+01 2e+01
## np 33.3% 33.3% 33.3%
## NAs 0 0 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 43.121, df = 2, p-value = 4.328e-10
## -------------------------------------------------------------------------
## ER_ugcm2h ~ treatment
##
## Summary:
## n pairs: 7e+01, valid: 7e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N P NP
## mean -2.3e+01 -2.1e+01 -2.5e+01 -2.4e+01
## median -2.5e+01 -1.7e+01 -2.0e+01 -1.8e+01
## sd 2.0e+01 1.1e+01 1.4e+01 1.4e+01
## IQR 3.8e+01 2.0e+01 2.4e+01 2.4e+01
## n 2e+01 2e+01 2e+01 2e+01
## np 25.0% 25.0% 25.0% 25.0%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 0.83723, df = 3, p-value = 0.8405
From our preliminary glance at the data, it looks like ER might be higher on average in the warmest temperature treatment (16 degrees). As for the nutrients, the control has the widest spread in ER, but the highest mean ER is in the N treatment.
## Anova Table (Type III tests)
##
## Response: ER_ugcm2h
## Sum Sq Df F value Pr(>F)
## temperature 14816.8 3 112.4701 < 2.2e-16 ***
## treatment 1500.2 3 11.3876 5.202e-06 ***
## temperature:treatment 2357.7 6 8.9481 5.598e-07 ***
## Residuals 2634.8 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), we see there is a evidence of a significant interaction effect (F(6,60)=8.9481, p<0.0001). The tests for the main effects of temperature treatment is significant (F(2,60)=63.4, p<0.00001) and so is nutrient treatment (F(3,60)=11.3876, p<0.0001.)
After fitting an ANOVA model it is important to always check the relevant model assumptions. This includes making QQ-plots and residual plots. Will use Levene’s test of equality of variances too.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 11 1.8059 0.07288 .
## 60
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Looks like we meet the assumptions of normality so we can move onto pairwise comparisons.
Pairwise comparsions for each level of temperature. At 8 degrees, there is a difference between the control and all nutrient treatments, but not between nutrient treatments. At 12 degrees, there is a different between C - N and C - NP. At 16 degrees, there is no statistically significant relationship between nutrient treatments.
## $emmeans
## temperature = 8:
## treatment emmean SE df lower.CL upper.CL
## C 0.5166667 2.705348 60 -4.894835 5.928168
## N -13.5333333 2.705348 60 -18.944835 -8.121832
## P -13.9833333 2.705348 60 -19.394835 -8.571832
## NP -21.3333333 2.705348 60 -26.744835 -15.921832
##
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C -26.4500000 2.705348 60 -31.861502 -21.038498
## N -14.3333333 2.705348 60 -19.744835 -8.921832
## P -18.7833333 2.705348 60 -24.194835 -13.371832
## NP -10.7833333 2.705348 60 -16.194835 -5.371832
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C -42.0666667 2.705348 60 -47.478168 -36.655165
## N -35.2166667 2.705348 60 -40.628168 -29.805165
## P -42.5833333 2.705348 60 -47.994835 -37.171832
## NP -39.7833333 2.705348 60 -45.194835 -34.371832
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 8:
## contrast estimate SE df t.ratio p.value
## C - N 14.0500000 3.82594 60 3.672 0.0028
## C - P 14.5000000 3.82594 60 3.790 0.0019
## C - NP 21.8500000 3.82594 60 5.711 <.0001
## N - P 0.4500000 3.82594 60 0.118 0.9994
## N - NP 7.8000000 3.82594 60 2.039 0.1855
## P - NP 7.3500000 3.82594 60 1.921 0.2302
##
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N -12.1166667 3.82594 60 -3.167 0.0126
## C - P -7.6666667 3.82594 60 -2.004 0.1980
## C - NP -15.6666667 3.82594 60 -4.095 0.0007
## N - P 4.4500000 3.82594 60 1.163 0.6522
## N - NP -3.5500000 3.82594 60 -0.928 0.7900
## P - NP -8.0000000 3.82594 60 -2.091 0.1678
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N -6.8500000 3.82594 60 -1.790 0.2880
## C - P 0.5166667 3.82594 60 0.135 0.9991
## C - NP -2.2833333 3.82594 60 -0.597 0.9327
## N - P 7.3666667 3.82594 60 1.925 0.2285
## N - NP 4.5666667 3.82594 60 1.194 0.6332
## P - NP -2.8000000 3.82594 60 -0.732 0.8839
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Pairwise comparsions for each level of nutrient treatment. In the control treatment, there is a difference between 8 - 12, 8 - 16, and between 12 - 16. In the N treatment, there is no difference between 8 and 12, but there is a difference between 8 - 16 and 12 - 16. In the P treatment, there is a difference between 8 - 16 and 12- 16. Lastly, in the NP treatment, there is a difference between 8 - 12, 8 - 16, and 12 - 16.
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 8 0.5166667 2.705348 60 -4.894835 5.928168
## 12 -26.4500000 2.705348 60 -31.861502 -21.038498
## 16 -42.0666667 2.705348 60 -47.478168 -36.655165
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 8 -13.5333333 2.705348 60 -18.944835 -8.121832
## 12 -14.3333333 2.705348 60 -19.744835 -8.921832
## 16 -35.2166667 2.705348 60 -40.628168 -29.805165
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 8 -13.9833333 2.705348 60 -19.394835 -8.571832
## 12 -18.7833333 2.705348 60 -24.194835 -13.371832
## 16 -42.5833333 2.705348 60 -47.994835 -37.171832
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 8 -21.3333333 2.705348 60 -26.744835 -15.921832
## 12 -10.7833333 2.705348 60 -16.194835 -5.371832
## 16 -39.7833333 2.705348 60 -45.194835 -34.371832
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 8 - 12 26.96667 3.82594 60 7.048 <.0001
## 8 - 16 42.58333 3.82594 60 11.130 <.0001
## 12 - 16 15.61667 3.82594 60 4.082 0.0004
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 8 - 12 0.80000 3.82594 60 0.209 0.9762
## 8 - 16 21.68333 3.82594 60 5.667 <.0001
## 12 - 16 20.88333 3.82594 60 5.458 <.0001
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 8 - 12 4.80000 3.82594 60 1.255 0.4263
## 8 - 16 28.60000 3.82594 60 7.475 <.0001
## 12 - 16 23.80000 3.82594 60 6.221 <.0001
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 8 - 12 -10.55000 3.82594 60 -2.757 0.0208
## 8 - 16 18.45000 3.82594 60 4.822 <.0001
## 12 - 16 29.00000 3.82594 60 7.580 <.0001
##
## P value adjustment: tukey method for comparing a family of 3 estimates
| pigment_ID | treatment | temperature | mean_ug_cm2 | sd_mass_ug_cm2 | n |
|---|---|---|---|---|---|
| chla | C | 12 | 2.422 | 1.358 | 3 |
| chla | N | 12 | 1.370 | 1.690 | 3 |
| chla | NP | 12 | 1.425 | 0.572 | 3 |
| chla | P | 12 | 2.688 | 1.682 | 3 |
| chla | C | 16 | 2.440 | 0.744 | 3 |
| chla | N | 16 | 2.866 | 0.815 | 3 |
| chla | NP | 16 | 1.926 | 0.662 | 3 |
| chla | P | 16 | 1.801 | 1.386 | 3 |
## -------------------------------------------------------------------------
## chla ~ temperature
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 2
##
##
## 12 16
## mean 2.0e+00 2.3e+00
## median 1.7e+00 2.3e+00
## sd 1.3e+00 9.2e-01
## IQR 1.5e+00 1.1e+00
## n 1e+01 1e+01
## np 50.0% 50.0%
## NAs 0 0
## 0s 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 0.75, df = 1, p-value = 0.3865
## -------------------------------------------------------------------------
## chla ~ treatment
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N NP P
## mean 2.4e+00 2.1e+00 1.7e+00 2.2e+00
## median 2.2e+00 2.4e+00 1.6e+00 1.9e+00
## sd 9.8e-01 1.4e+00 6.2e-01 1.5e+00
## IQR 1.1e+00 2.1e+00 6.6e-01 1.9e+00
## n 6e+00 6e+00 6e+00 6e+00
## np 25.0% 25.0% 25.0% 25.0%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 1.0133, df = 3, p-value = 0.798
## Anova Table (Type III tests)
##
## Response: chla
## Sum Sq Df F value Pr(>F)
## (Intercept) 17.6031 1 12.3167 0.002905 **
## temperature 0.0004 1 0.0003 0.986164
## treatment 4.1320 3 0.9637 0.433899
## temperature:treatment 4.4382 3 1.0351 0.403617
## Residuals 22.8673 16
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), <
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 7 0.3228 0.9325
## 16
Looks like we meet the assumptions of normality, but no need to move onto pairwise comparisons.
| pigment_ID | treatment | temperature | mean_ug_cm2 | sd_mass_ug_cm2 | n |
|---|---|---|---|---|---|
| chlb | C | 12 | 0.933 | 0.488 | 3 |
| chlb | N | 12 | 0.573 | 0.654 | 3 |
| chlb | NP | 12 | 0.656 | 0.294 | 3 |
| chlb | P | 12 | 1.074 | 0.588 | 3 |
| chlb | C | 16 | 1.004 | 0.244 | 3 |
| chlb | N | 16 | 1.007 | 0.148 | 3 |
| chlb | NP | 16 | 0.912 | 0.225 | 3 |
| chlb | P | 16 | 0.499 | 0.248 | 3 |
## -------------------------------------------------------------------------
## chlb ~ temperature
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 2
##
##
## 12 16
## mean 8.1e-01 8.6e-01
## median 7.6e-01 8.9e-01
## sd 4.9e-01 2.9e-01
## IQR 5.7e-01 2.6e-01
## n 1e+01 1e+01
## np 50.0% 50.0%
## NAs 0 0
## 0s 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 0.40333, df = 1, p-value = 0.5254
## -------------------------------------------------------------------------
## chlb ~ treatment
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N NP P
## mean 9.7e-01 7.9e-01 7.8e-01 7.9e-01
## median 8.6e-01 9.2e-01 8.1e-01 6.9e-01
## sd 3.5e-01 4.9e-01 2.7e-01 5.1e-01
## IQR 3.4e-01 6.8e-01 2.9e-01 3.9e-01
## n 6e+00 6e+00 6e+00 6e+00
## np 25.0% 25.0% 25.0% 25.0%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 0.9, df = 3, p-value = 0.8254
## Anova Table (Type III tests)
##
## Response: chlb
## Sum Sq Df F value Pr(>F)
## (Intercept) 2.61308 1 16.1806 0.0009843 ***
## temperature 0.00753 1 0.0466 0.8317406
## treatment 0.49485 3 1.0214 0.4092645
## temperature:treatment 0.87158 3 1.7990 0.1879399
## Residuals 2.58391 16
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), <
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 7 0.4584 0.8503
## 16
Looks like we meet the assumptions of normality, but no need to move onto pairwise comparisons.
| pigment_ID | treatment | temperature | mean_ug_cm2 | sd_mass_ug_cm2 | n |
|---|---|---|---|---|---|
| perc_green | C | 12 | 39.467 | 3.370 | 3 |
| perc_green | N | 12 | 45.440 | 7.199 | 3 |
| perc_green | NP | 12 | 45.235 | 6.726 | 3 |
| perc_green | P | 12 | 42.339 | 7.957 | 3 |
| perc_green | C | 16 | 41.803 | 4.855 | 3 |
| perc_green | N | 16 | 35.982 | 4.384 | 3 |
| perc_green | NP | 16 | 48.411 | 6.013 | 3 |
| perc_green | P | 16 | 31.972 | 8.403 | 3 |
## -------------------------------------------------------------------------
## perc_green ~ temperature
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 2
##
##
## 12 16
## mean 4.3e+01 4.0e+01
## median 4.1e+01 3.9e+01
## sd 6.1e+00 8.3e+00
## IQR 9.1e+00 6.2e+00
## n 1e+01 1e+01
## np 50.0% 50.0%
## NAs 0 0
## 0s 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 0.56333, df = 1, p-value = 0.4529
## -------------------------------------------------------------------------
## perc_green ~ treatment
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N NP P
## mean 4.1e+01 4.1e+01 4.7e+01 3.7e+01
## median 3.9e+01 3.9e+01 4.8e+01 3.8e+01
## sd 4.0e+00 7.4e+00 6.0e+00 9.3e+00
## IQR 4.2e+00 3.8e+00 8.3e+00 5.1e+00
## n 6e+00 6e+00 6e+00 6e+00
## np 25.0% 25.0% 25.0% 25.0%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 5.4667, df = 3, p-value = 0.1406
## Anova Table (Type III tests)
##
## Response: perc_green
## Sum Sq Df F value Pr(>F)
## (Intercept) 4672.9 1 116.3567 9.48e-09 ***
## temperature 8.2 1 0.2038 0.6577
## treatment 71.4 3 0.5928 0.6287
## temperature:treatment 241.9 3 2.0076 0.1535
## Residuals 642.6 16
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Studying the output of the ANOVA table (Type III SS), <
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 7 0.1917 0.983
## 16
Looks like we meet the assumptions of normality, but no need to move onto pairwise comparisons.
## Observations: 72
## Variables: 19
## $ treatment <chr> "C", "C", "C", "C", "C", "C", "N", "N", "N", "N"...
## $ temperature <chr> "8", "8", "8", "8", "8", "8", "8", "8", "8", "8"...
## $ sample_id <chr> "8C-1", "8C-2", "8C-3", "8C-4", "8C-5", "8C-6", ...
## $ rep <int> 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, ...
## $ NEP_ugcm2h <dbl> 65.9, 28.8, 24.1, 26.6, 26.2, 54.0, 31.3, 29.4, ...
## $ ER_ugcm2h <dbl> 0.1, 6.5, 5.2, -0.9, -4.0, -3.8, -15.5, -9.7, -1...
## $ GPP_ugcm2h <dbl> 66.1, 35.3, 29.3, 27.4, 30.1, 57.8, 46.8, 39.0, ...
## $ changePO4_uM <dbl> -6.681, -6.078, -5.923, -6.597, -1.115, -6.336, ...
## $ changeTDN_uM <dbl> 5.887, 5.967, 9.587, 10.057, 6.307, 5.417, 20.59...
## $ changeDOC_mgL <dbl> 3.281, 3.241, 2.091, 2.781, 3.271, 1.891, 3.001,...
## $ changeNH4_mgL <dbl> -0.032, -0.022, 0.018, 0.018, -0.012, -0.012, -0...
## $ percTDNuptake <dbl> 32.54824, 32.99055, 53.00492, 55.60347, 34.87035...
## $ percPO4uptake <dbl> -835.125000, -759.750000, -740.375000, -824.6250...
## $ chl_a_phyto <dbl> 23.52, 11.44, 23.30, 25.36, 22.91, 13.15, 14.74,...
## $ DOC_mgL <dbl> 3.94, 3.90, 2.75, 3.44, 3.93, 2.55, 3.66, 2.58, ...
## $ NH4_mgL <dbl> 0.07, 0.06, 0.02, 0.02, 0.05, 0.05, 0.26, 0.31, ...
## $ chla <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, ...
## $ chlb <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, ...
## $ perc_green <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, ...
## cor
## 0.3608053
## cor
## 0.3410461
## cor
## -0.4436344
## cor
## 0.3945505
## cor
## 0.3761142
## cor
## -0.4659032
## cor
## 0.8316538
## cor
## 0.2653822
## cor
## -0.619436
## cor
## -0.4367146
## cor
## -0.3221464
## cor
## 0.256881
## -------------------------------------------------------------------------
## BPratio ~ temperature
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 2
##
##
## 12 16
## mean 7.2e-02 1.0e-01
## median 4.4e-02 7.1e-02
## sd 5.8e-02 7.5e-02
## IQR 8.1e-02 9.8e-02
## n 1e+01 1e+01
## np 50.0% 50.0%
## NAs 0 0
## 0s 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 1.92, df = 1, p-value = 0.1659
## -------------------------------------------------------------------------
## BPratio ~ treatment
##
## Summary:
## n pairs: 2e+01, valid: 2e+01 (100.0%), missings: 0 (0.0%), groups: 4
##
##
## C N NP P
## mean 1.4e-01 6.6e-02 3.0e-02 1.2e-01
## median 1.2e-01 6.3e-02 2.8e-02 1.0e-01
## sd 7.4e-02 3.9e-02 1.3e-02 7.2e-02
## IQR 8.6e-02 4.7e-02 1.5e-02 1.1e-01
## n 6e+00 6e+00 6e+00 6e+00
## np 25.0% 25.0% 25.0% 25.0%
## NAs 0 0 0 0
## 0s 0 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 11.107, df = 3, p-value = 0.01116
Lower B:P means less benthic biomass relative to planktonic biomass
## Anova Table (Type III tests)
##
## Response: BPratio
## Sum Sq Df F value Pr(>F)
## (Intercept) 0.025509 1 8.6099 0.009724 **
## temperature 0.012429 1 4.1952 0.057317 .
## treatment 0.018302 3 2.0591 0.146092
## temperature:treatment 0.008785 3 0.9884 0.423172
## Residuals 0.047404 16
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 7 0.6698 0.6948
## 16
Pairwise comparsions for each level of temperature. Only at 16 degrees C do we see a difference between C & NP.
## $emmeans
## temperature = 12:
## treatment emmean SE df lower.CL upper.CL
## C 0.09221148 0.03142578 16 0.02559180 0.15883117
## N 0.04924343 0.03142578 16 -0.01737625 0.11586311
## NP 0.02195727 0.03142578 16 -0.04466242 0.08857695
## P 0.12367432 0.03142578 16 0.05705464 0.19029401
##
## temperature = 16:
## treatment emmean SE df lower.CL upper.CL
## C 0.18323962 0.03142578 16 0.11661994 0.24985931
## N 0.08340603 0.03142578 16 0.01678634 0.15002571
## NP 0.03811816 0.03142578 16 -0.02850152 0.10473785
## P 0.10966228 0.03142578 16 0.04304260 0.17628196
##
## Confidence level used: 0.95
##
## $contrasts
## temperature = 12:
## contrast estimate SE df t.ratio p.value
## C - N 0.04296805 0.04444277 16 0.967 0.7698
## C - NP 0.07025421 0.04444277 16 1.581 0.4164
## C - P -0.03146284 0.04444277 16 -0.708 0.8925
## N - NP 0.02728616 0.04444277 16 0.614 0.9261
## N - P -0.07443089 0.04444277 16 -1.675 0.3681
## NP - P -0.10171706 0.04444277 16 -2.289 0.1423
##
## temperature = 16:
## contrast estimate SE df t.ratio p.value
## C - N 0.09983359 0.04444277 16 2.246 0.1530
## C - NP 0.14512146 0.04444277 16 3.265 0.0226
## C - P 0.07357734 0.04444277 16 1.656 0.3777
## N - NP 0.04528787 0.04444277 16 1.019 0.7410
## N - P -0.02625625 0.04444277 16 -0.591 0.9334
## NP - P -0.07154412 0.04444277 16 -1.610 0.4011
##
## P value adjustment: tukey method for comparing a family of 4 estimates
## $emmeans
## treatment = C:
## temperature emmean SE df lower.CL upper.CL
## 12 0.09221148 0.03142578 16 0.02559180 0.15883117
## 16 0.18323962 0.03142578 16 0.11661994 0.24985931
##
## treatment = N:
## temperature emmean SE df lower.CL upper.CL
## 12 0.04924343 0.03142578 16 -0.01737625 0.11586311
## 16 0.08340603 0.03142578 16 0.01678634 0.15002571
##
## treatment = NP:
## temperature emmean SE df lower.CL upper.CL
## 12 0.02195727 0.03142578 16 -0.04466242 0.08857695
## 16 0.03811816 0.03142578 16 -0.02850152 0.10473785
##
## treatment = P:
## temperature emmean SE df lower.CL upper.CL
## 12 0.12367432 0.03142578 16 0.05705464 0.19029401
## 16 0.10966228 0.03142578 16 0.04304260 0.17628196
##
## Confidence level used: 0.95
##
## $contrasts
## treatment = C:
## contrast estimate SE df t.ratio p.value
## 12 - 16 -0.09102814 0.04444277 16 -2.048 0.0573
##
## treatment = N:
## contrast estimate SE df t.ratio p.value
## 12 - 16 -0.03416260 0.04444277 16 -0.769 0.4533
##
## treatment = NP:
## contrast estimate SE df t.ratio p.value
## 12 - 16 -0.01616089 0.04444277 16 -0.364 0.7209
##
## treatment = P:
## contrast estimate SE df t.ratio p.value
## 12 - 16 0.01401204 0.04444277 16 0.315 0.7566